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Abstract 

A new finite element for coupling built-up shell 
substructures is presented. The present work extends 
the hybrid variational formulation of the interface ele- 
ment developed by Aminpour and Ransom to permit 
coupling between two intersecting substructures. 
Designed for the assembly of independent built-up finite 
element models, this technique provides a level of mod- 
eling flexibility previously unavailable. 

Introduction 

The analysis of modern aerospace and ground 
vehicle structures relies heavily on the finite element 
method. While detailed and accurate stress predictions 
can be made for even the most geometrically complex 
structural assemblies, practical limitations in mesh size 
and mesh discretization often prevent the method’s use 
in the design phase. Indeed, the finite element method 
can overwhelm available resources in two major ways. 
Not only can the size of a practical model quickly reach 
the limits of the accessible computational hardware, but 
the generation and validation of the mesh are likely to 
require a more than considerable effort. In practice, 
many useful analysis and design studies that could lead 
to lower-cost designs are precluded. 

One example where meshing can become extraordi- 
narily costly is structural analysis in multi-disciplinary 
optimization. Shape optimization in aircraft design 
involves iterations on aerodynamic and structural 
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analyses. These optimizations are carried out concur- 
rently, with both analyses sharing a common descrip- 
tion of the airframe. This parameterized three-dimen- 
sional description is used to develop the aero grid and 
the structural finite element mesh of the skin. 
Although similar to each other, these two meshes are so 
fundamentally different in purpose that the same dis- 
cretization is not used for aerodynamic and structural 
analyses. Indeed, aero grids are only discretizations of 
the boundary of the aircraft's volume, and their size and 
density is often changed during an analysis. The struc- 
tural mesh of the skin, on the other hand, is tightly 
constrained by the locations of its internal components. 
Node-to-node compatibility between structural compo- 
nents such as bulkheads, ribs, stiffeners and the aircraft 
skin dictate constraints on the finite element mesh that 
do not have a counterpart in the aero grid. 
Consequently, structural skin meshes are, at best, diffi- 
cult and expensive to generate. 


Figure 1. The generation of structural mod- 
els would be streamlined if the wing box and 
ribs could be connected to the aero grid. 

On the other hand, if the internal structural compo- 
nents of a wing could be coupled to the elements of the 
skin along arbitrary intersection curves, as illustrated in 




Fig. 1, the constraints that make structural meshes dif- 
ficult to generate would be eliminated. Furthermore, 
the aero grid may now be suitable as a structural mesh 
for the skin. This would eliminate the need for remesh- 
ing, and provide the added advantage of relieving the 
translation losses that result from interpolating the 
aerodynamic loads from the aero grid to the structural 
mesh. 

Detailed finite element models are not efficient in 
the design phase because they are computationally 
intensive and require extensive modeling effort. The 
interface element developed by Aminpour et al. 1 and 
Ransom et al. 2 addresses some of the difficulties 
associated with mesh generation and optimization. This 
interface element was developed to enforce compatibility 
on the boundaries of independently modeled shell 
substructures so that the need for complex mesh 
transitioning is eliminated. Because it allows the finite 
element domain to be subdivided into various 
substructures with different element discretizations, the 
interface element provides modeling flexibility and 
allows the generation of independently optimal meshes. 

Unlike other coupling methods based on mortar el- 
ements 3 or spline fit 4 , the interface element based on a 
hybrid variational formulation has been shown to be ro- 
bust and accurate. It preserves solution accuracy of dis- 
placements and stresses across the common substructure 
boundaries. 

However, the interface element of Refs. 1 and 2 
lacks generality in one essential respect. It is designed 
for coupling shell substructures sharing a common edge 
on which all the connected finite element nodes of the 
connected substructures lie (Figs. 2a and 2b). The pre- 
sent formulation extends this concept to include the 
more general problem of coupling across element faces, 
as shown in Fig 3. In this figure, substructure £2] is 
coupled across the surface of substructure £2 2 . For sim- 
plicity, we refer hereafter to £2j as the stiffener, and to 
£2 2 as the skin. 



Figure 2. The interface element of Refs. 1 
and 2 allows coupling of independently mod- 
eled substructures along a common edge. 



contributing to the interlace 

Figure 3. A new cross-surface interface 
element is used to couple substructure £2| 
across the surface of substructure £2 2 . 

The objective of this work is to develop an inter- 
face element with cross-surface capabilities, and to 
evaluate the accuracy and usefulness of this element. 
Therefore, it was considered of vital importance that 
even the prototype for this element possess a fully au- 
tomated capability for determining the path of the 
stiffener on the surface of the skin. This capability in- 
cludes the determination of the elements and degrees of 
freedom that contribute to the interface. The first two 
examples presented in this paper illustrate the accuracy 
of the method under a ’’worst-case scenario”, as well as 
under more normal conditions. The last problem is an 
example of a parametric study in which the location of 
structural components is varied without modifying the 
mesh. 

Cross-Surface Interface Element 

The hybrid variational formulation used in the in- 
terface element was presented in detail by Aminpour et 
al. 1 and is briefly described here. The formulation em- 
ploys an integral form for the compatibility between the 
interface and the subdomains. For two subdomains £2j 
and £2 2 sharing a common interface, 5, the modified to- 
tal potential energy is written as 

n = n Q +n n \](v-u,)+?J 2 (v-u t )ds (i) 

where Aj and are vectors of the Lagrange multipliers 
corresponding to subdomains £2] and £2 2 , respectively. 
The displacement vector of the interface element, v , is 
allowed to be independent of the displacement vectors 
Uj and « 2 , of the subdomains. It can be shown that the 
stationary condition for the total potential energy results 
in the Euler equations 

[ = (<*«), . '='< 2 

Aj + A 2 = 0 

", = v 
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sn = o => 
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on s (2) 


where n represents the outward unit normal to the sub- 
domains along their interfaces, a is the stress tensor, 
and the Lagrange multipliers Aj and represent the 
tractions at the interfaces. Equation 2 is an expression 
of displacement continuity and traction reciprocity, 
which guarantees that equilibrium is maintained across 
the interface. 

The interpolations for the finite element displace- 
ments, interface displacements, and interface tractions 
are, respectively, 

= N i <ii 

V = T q s * i = 7, 2 for Q, and ft 2 (3) 

A, = /?, A, . 

where q { , and q s are the nodal degrees of freedom 
corresponding to u i andv, respectively, are the 
unknown coefficients for the Lagrange multipli- 
ers, , Ni is the matrix of the finite element shape 
functions on substructure Q n T is the matrix of shape 
functions for the (independent) interface displacements, 
and Rj is the matrix of interpolating functions for the 
Lagrange multipliers on substructure Q r 

Taking the first variation of the modified total po- 
tential energy, Eq. 1, with respect to the independent 
degrees of freedom, and setting the variation to zero 
yields the system of equations 

~K X 0 0 A/, Olfftj f// 

0 K 2 0 0 M 2 q 2 f 2 

0 0 0 Gj G 2 < q s ► = < 0 " (4) 

M t x 0 G? 0 0 A, 0 

0 M t 2 Gl 0 oJ[A 2 J [o 

where Kj and K 2 are the assembled stiffness matrices of 
the uncoupled subdomains, are the nodal force vec- 
tors, and the coupling submatrices M ( and G t are com- 
posed of integrals of the interpolation matrices N r T, 
and R { along the interface, 5. They are defined as 

M, = -j>, r R t ds and G, = [t j R, ds (5) 

For simplicity, the following development is for a 
single degree-of-freedom. The same pattern is repeated 
for the remaining degrees of freedom. Furthermore, the 
following computation of the integrals in Eqs. 5 is 
shown for one four-node skin element along the path of 
the stiffener. The full matrices are obtained by adding 
the terms for all the elements crossed by the stiffener. 
The computation of the coupling terms for elements 
other than four-node is not attempted here. 


For a four-node element, n 2 , n v n 4 }., 

where the bilinear interpolation functions are 

« 2=±0 + C ) 0 - r ?) (6) 

n, = i(l + £)0 + *0 

The distribution of the tractions, X t , is assumed to 
be constant over the face of a four-node element, so the 
corresponding interpolation function matrix is /?• = 1. 
Therefore, the form of the matrix M , is 
A/. = {m p m 2 , m 3 , m 4 }. , where m } is a function of n h 

For the original interface element, illustrated in 
Fig. 4a, the integration of matrix A// in Eq. 5 is only 
in rj, from rj=-l to rj=+7 . Substituting £=-/ into 
Eqs. 6, and noting that ds = 1 drf, the terms in the 

coupling matrix M t are 

(7) 

m 3 = — J- J J *0-(1 + tj) - 1 ds = 0 

Jj\ 

m 4 =-ip(l + ri)lds = -&^ 

,/s i z 

which confirms that the degrees of freedom at nodes 2 
and 3 (corresponding to m 2 and m 3 ) do not contribute to 
the interface. 

In the general case shown in Fig. 4b, however, the 
interface crosses the finite element's face, and all nodes 
in the element contribute to the interface. The interpo- 
lation is expressed in terms of both natural coordinates, 
rj and £ The first term of M x is 


m, = -i p (i - cxi - n) - 1 ^ w 

Js\ 



Figure 4a. Integration path along the edge of 
an element. 
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finite element nodes 



Figure 4b* Integration path across the sur- 
face of an element. 

To perform the integration in Eq. 8 it is necessary 
to relate the element's natural coordinates, f and 77 , and 
the path increment, ds, to a single variable. This 
requires a complete description of the intersection path 
that is generated by the stiffener on the elements in the 
skin. One simplifying assumption that allows closed- 
form solution of this integration is illustrated in Fig. 5. 
It is assumed that the interface is contained in a plane 
which is referred to as the plane of the interface. 
The equation of this plane is expressed in terms of its 
normal, n-{n x ,ny,n z } and offset, d 

n x x + n y y + n.z = d (9) 



Figure 5. A closed-form solution is ob- 
tained by assuming that the interface is con- 
tained in a plane which intersects the skin 
(substructure 0 2 ). 

The first step in the automated definition phase of 
the interface element consists of computing the four co- 
efficients in Eq. 9 from the nodes on the stiffener 


(provided by the user) which define the location of the 
interface. The interface plane is used to locate the ele- 
ments in the skin that lie along the path of the stiffener. 
For these elements, two and only two of the following 
four intersection conditions are valid: 

77 = — 1 , -1 < £< 1 

C=+i, - 1 < 77 < 1 
rj = +l, -1 < £<1 

f =-l> — 1 < 77 < 1 

For instance, on the edge between nodes 1 and 2 in Fig. 
4b, 77=— 1 (first condition in Eqs. 10). Substituting 
this condition into the interpolation functions, Eqs. 6 , 
and then into the equation of the intersecting plane, Eq. 
9, and solving for £ results in 

^_ ^(- t l +^) + » v (>’ l +> 2 ) + ^(^ 1 + z 2 )- 2 </ 

b «x(^-^) + nvU->’2) + « z (zi-Z 2 ) ' 

If -1<£<+1, the first condition in Eqs. 10 is 
satisfied. One more condition out of the remaining 
three must also be satisfied for a valid intersection. 
After sorting the elements along the intersection of the 
plane and the skin, it is possible to identify and 
eliminate the elements outside the path of the interface. 
All of the remaining elements contribute to the 
interface, and they possess a point through which the 
path of the interface enters the element and a point 
through which it exits the element. We refer to the 
natural coordinates of these two points as Q ,r\ l and 
T)°, respectively, and to the corresponding Cartesian 
coordinates as p‘ = { p ‘ x , p ‘ y , p [ } and p° = f p ° , p ° v , p " } 
(see Fig. 4b). 

The relation between £ and 77 needed to perform the 
integration of Eq. 8 is obtained as follows. The coordi- 
nates of a point on an element can be written using the 
element's interpolation as: 

4 x = £ri(x ] -x 2 +Xi-x 4 )+£(-x } +x 2 + * 3 -jr 4 ) 

+ t](-x l - x 2 + JC 3 + x 4 ) + (jCj + x 2 + *3 -I- x 4 ) 

4 >' = ^(>’1 - >’ 2 + y, - >4 ) + C(-.v, + >’ 2 + y 3 -y 4 ) 

+ t](-y, - >’ 2 + y 3 + y A ) + (>’i + y 2 + y 3 + y 4 ) 

4z = frj(z, - z 2 + z, - z 4 ) + £(-Z| +z 2 + z, - z 4 ) 

+ »7(-Z| -z 2 + z 3 + z 4 ) + (z, + z 2 + z, + z 4 ) 


( 10 ) 


where , y, , and z n i = 1,4, are the coordinates of the 
nodes of the element. Substituting Eq. 12 into the 
equation of the plane of the interface, Eq. 9, gives 



where the constants r, are defined: 

t, = « t (*i -x 2 +x^-x 4 ) + n v (y, - y 2 + y, - y 4 ) 

+ n,(z t -z 2 + z, -z 4 ) 

U = n x (-jr, + x 2 + x, - x 4 ) + n y (— y, + y 2 + y, - y 4 ) 

+ n.(-z, + z 2 + z, -z 4 ) 

' (14) 

h = ", (~ x \ ~x 2 +x,+x i ) + n y (-y, - y 2 + y } + y 4 ) 

+ rt.(-z, -z 2 + z, + z 4 ) 

f 4 = n x ( ■ x I + x 2 + *1 + x 4 ) + n v ( y, + y 2 + y 3 + y 4 ) 

+ n. (z, + z 2 + z, + z 4 ) - 4 d 

Finally, the incremental variable ds in Eq. 8 can be 
expressed in terms of the incremental coordinates as 

dx = xjdx 2 + dy 2 + dz" (15) 

In the case of a four-node element, we can write the 
linear relations between the Cartesian intersection coor- 
dinates and the path increments as 

dx _ dy _ dz 

(k-k)"(k-k)"(p=-p;) 

from which two of the increments may be expressed in 
terms of a third. For instance, if (p° x - p‘ x ) in Eq. 16 is 
nonzero, then Eq. 15 gives 

Js = Vl + a 2 +)3 2 \dx\ (17) 


where A, B, C, D and E are constants which depend 
only on the element's nodal coordinates. Similar results 
are obtained for the remaining three terms, m 2 . m 3 , and 
m 4 . 

Results and Discussion 

The following examples were developed to evaluate 
the accuracy and limitations of the cross-surface inter- 
face element. A benchmark problem is presented first 
to assess the performance of the elements under a "worst 
case" condition, i.e., when the interface produces a se- 
vere discontinuity in the shell elements under the path 
of the interface. The second example illustrates the si- 
multaneous use of four interface elements to couple 
blade stiffeners to a composite flat plate. Finally, the 
last problem consists of a framed fuselage panel with 
two curved interfaces. The use of interface elements in 
parametric studies is illustrated with two solutions cor- 
responding to different placements of the frames. 

Case 1 : Cantilevered plate with stiffener supp ort 

The most demanding problems for cross-surface 
coupling involve a severe discontinuity in the transverse 
forces of the shell elements in the skin. The following 
"worst-case scenario" was designed to introduce a 
discontinuity which the shell elements in the skin 
cannot accomodate. 


where a and (5 are the appropriate ratios of the denomi- 
nators in Eqs. 16. The term dx is obtained from the 
First term in Eqs. 12. Therefore, we can write 


A = VEZig 

4 


( 


O] -x 2 ) 


(1 TJ )( f 2 6 f 4 ) 


“(*3 --* 4 ) 

“(*1 -* 4 ) 
“(*2 -* 3 ) 


(l + 77)(r 2 f,-f,r 4 

(^i+'2) 2 


l n*i+t 2 ) 

\ \ 

i_l£i±k 
y *7 U + h j j 


) 


dr] 


(18) 


Using (18), and (13), it can be shown that the closed 
form solution for the integral in Eq. 8 is 


m | = A(1 1,-TJi) + B(T]l-T]f) 


+ C 


(/ 2 +r,?7,.) 2 (f 2 +r,J7„) 


(19) 


+ D 


1 


1 


\(t 2 +'.*?, ) (fj+W.) 


+ E log 


h±Lik 

h+wj 


This benchmark problem consists of a cantilevered 
narrow plate. The plate is supported at 3/4 of its length 
by a transverse stiffener which is clamped at its free end, 
as shown in Fig. 6. A load of 1.0 lb / in is applied at one 
end of the plate. Young's modulus is E=10 6 psi., the 
thickness is 0.1 in., and Poisson's ratio is taken to be 
v=0 to achieve pure cylindrical bending. 

The mesh of the stiffener consists of three four- 
node elements along the width and two along the 
height. When the number of elements along the length 
of the plate is not a multiple of four, the elements under 
the stiffener are in the region of the model where the 
bending moment changes from its maximum value on 
the left of the stiffener, to zero, on the right of the stiff- 
ener. Since the four-node elements used here can only 
represent a constant bending moment M x and transverse 
shear force Q x along the x axis, approximation errors re- 
sult, as shown in Fig. 6. This case represents the most 
stringent interface problem, since the stiffener 
introduces a severe discontinuity in the moment and 
transverse force for the elements on the plate. 
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The plate's tip deflection is shown in Fig. 7 versus 
the number of elements along the span of the plate. 
When the number of elements is a multiple of four, the 
stiffener intersects along the edges of elements in the 
plate, and the solution is equivalent to the reference so- 
lution (i.e., model with no interface element). 
However, the convergence is not monotonic, and when 
the stiffener intersects near the center of an element on 
the plate, as for N=6, the error reaches a maximum 
(12%). It can also be observed that the solution is al- 
ways bracketed between the two closest standard finite 
element models (i.e., the models with no interface ele- 
ment). The lower bracket consists of a model where the 
stiffener is relocated to the left edge of the intersected 
shell elements, and the upper bracket is one where the 
stiffener is placed at the right edge of the intersected 
shell element. As expected, the bracket narrows as the 
mesh in the plate is refined. It can be observed that the 
error in the coupled model is always much smaller than 
the bracket. In other words, the error introduced by the 
interface element is much smaller than the error that 
would result from attaching the stiffener to the closest 
element edge. 

It is important to note that the oscillation in the er- 
ror is not caused by the interface element, but rather, by 
the finite element discretization of the plate. A similar 
error would occur, for instance, if one were to place a 


line load across the shell elements in a standard finite 
element mesh with no interface. 

tip deflection 
(inches) 



number of elements along span 

Figure 7. Tip displacements for cantilevered 
plate with stiffener support. The convergence of 
the solution is not monotonic. 
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Case 2: Blade-stiffened composite panel with a hole 

A linear stress analysis of a flat, blade-stiffened 
graphite-epoxy panel loaded in axial compression is per- 
formed. The panel is 30 in. long, the width is 1 1.5 in., 
the stiffener spacing is 4.5 in., the stiffener height is 
1.4 in., and the cutout radius is 1 in. All three blade- 
stiffeners have identical cross-sections. The loading is 
uniform axial compression in the x direction. All 
degrees of freedom of the loaded ends of the panel are 
constrained except for the motion in the axial direction. 
The material system used for the panel is T300/5208 
graphite-epoxy unidirectional tapes with a nominal ply 
thickness of 0.0055 in. 

The panel skin is a 25-ply symmetric laminate 
([ ±45 / 0 2 / +45 / 0 3 / ±45 / 0 / 0]) s and the blade stiff- 
eners are 24-ply symmetric laminates ([±45/ 0, o ]) s . 

Since no theoretical solution exists for this con- 
figuration, a reference solution is obtained by using a 
finite element model which has no interface element. 
This reference finite element model is identical to that 
used to predict the global response of the panel in an 
earlier analysis by Ransom 5 . The finite element model 
used in the coupled analysis makes use of an annular 


mesh for the skin and rectangular meshes for the 
stiffeners. Interface elements are used to couple each of 
the stiffeners to the skin. The interface elements used 
in the outer two stiffeners intersect the faces of the fi- 
nite elements in the skin. 

Contour plots of the out-of-plane deflection, w, are 
shown in Fig. 8 superimposed on the deformed shapes 
of both the coupled model and the reference model. The 
distribution of the out-of-plane deflection at the panel 
mid-width is shown in Fig. 9 as a function of the longi- 
tudinal distance from one end of the panel normalized by 
the panel length. The distribution of the longitudinal 
stress resultant N x normalized by the applied load 
(N x ) aV g ( Le ., applied load divided by the panel width) is 
shown in Fig. 10 as a function of the lateral distance 
normalized by the cutout radius, R. Figures 9 and 10 
indicate that the coupled analysis is in good agreement 
with the reference solution. The coupled analysis yields 
slightly different values for the longitudinal stress 
resultants in the region around the hole and at the skin- 
stiffener intersection. These differences are partly due to 
the different discretizations used in the panel skin for the 
two analyses (The panel skin is slightly more refined in 
the coupled model than in the reference model). In any 
case, the results are in overall good agreement. 



a) Coupled Model 


b) Reference Model 


Max 


w 


Min 




Figure 8. Compression-loaded composite blade-stiffened panel: the out-of plane displacement con- 
tours indicate excellent agreement between the two models 
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Displacement, 

w 



thickness of the skin is 0.075 inches, and the following 
material properties were used: E=10 7 , v=0.3. The 
window in the panel is 1 1 inches in radius and is filled 
with a material that is 0.3 inches thick, with a Young’s 
modulus of E=5X10 6 , and a Poisson's ratio of v=0.3. 
The curved edges of the panel are clamped, while the 
straight edges and the free length of the frames are free 
in the axial (z) direction and constrained in the other 
directions. The finite element models for the frame 
were generated independently of that for the skin, and 
the frames may be placed at any position along the axis 
of the fuselage. Two frame spacings are presented here. 
In the first model, the frames are separated by a distance 
of 39.6 inches. In the second model, the spacing is 
increased to 48.6 inches. The resulting radial 
displacement contours are illustrated in Figs. 12a and 
12b. These models illustrate the potential use of the 
cross-surface interface elements in parametric studies of 
built-up structures. 


Figure 9. Normal displacement along the 
panel mid-width. 

N x 

Navg 



R 



Figure 11. Model of fuselage panel with a 
window and two frame. 


Figure 10. Axial stress distribution along 
the panel mid-length. 


Case 3: Fuselage panel with a window and two frames 

The final example shown in Fig. 1 1 consists of a 
fuselage panel with two circumferential frames. The 
diameter of the fuselage is 100 inches, the panel spans 
an arc of 34 degrees, and is 60 inches long. The 
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(a) (b) 

Figure 12. Radial displacement contours for two positions of the frames. 


Concluding Remarks 

This paper has presented the formulation of a new 
interface element with cross-surface capabilities. The 
objective of this work was to assess the potential of 
this technique in terms of its accuracy and, just as im- 
portantly, its ease of use. The implementation pre- 
sented here is only capable of coupling across four-node 
shell elements. These are the simplest to develop, yet 
their basic stress distributions make them also the least 
accurate. It is expected that elements of higher order 
will provide better accuracy by allowing a closer ap- 
proximation of the local discontinuities introduced by 
the interface. 

The results presented herein illustrate that the inter- 
face element can produce accurate results with even rela- 
tively coarse meshes. In fact, as the mesh is refined, 
any desired level of accuracy can be obtained. However, 
the results have also shown that the convergence is not 
always monotonic, and that the accuracy is reduced 
when the load path introduces a severe discontinuity in 
the transverse forces on the skin elements. 

A great deal of attention was given to developing 
this interface into a simple-to-use design tool for the 
quick generation of models. Therefore, even the proto- 
type implementation of this element includes the auto- 
matic computation of the intersection path of the stiff- 
ener across the skin. The user only specifies the por- 
tion of the stiffener along which the interface is to be 
placed (by specifying the sequence of node numbers on 
the stiffener edge). A closed-form solution for the path 
and coupling terms is obtained for a geometric restric- 
tion: it is assumed that the interface is contained in a 
plane that cuts through the skin. This limitation is not 
overly restrictive since substructures such as bulkheads, 
wing ribs, or blade stiffeners, are often plane. 
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